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Abstract. We summarize some recent results on the application of macroscopic 
spectral properties of random matrix models (RMM) to the QCD spectra. A 
comparison to existing lattice simulation is presented both for staggered and 
Wilson fermions for high but finite temperature. We consider two type of 
mixing between the four lowest Matsubara modes, corresponding to third and 
fifth order algebraic equation for the pertinent resolvent, respectively. 

1. Introduction 

Dedicated lattice simulations on the aspects of the QCD phase transition have 
triggered a number of theoretical investigations aimed at understanding some of the 
nonperturbative aspects of QCD whether in vacuum or matter Q. In particular, 
it was shown recently that the bulk character of the QCD Dirac spectra could be 
understood using simple random matrix models Q-Q, and that the spectral- level 
correlations are consistent with the general lore of random matrix theory 

At zero temperature but for massive quarks, QCD simulations with Wilson 
quarks show that the Dirac spectrum exhibits a structural change with increasing 
quark mass, a feature that is reflected in the random matrix version in terms of 
heavy quarks becoming localized over their Compton wavelength ||. A similar 
behavior is also observed at finite temperature when the quark modes are restricted 
to the lowest Matsubara modes ujq = = ttT J|,0]. 

At high temperature, the static screening lengths in QCD are dominated by 
the lowest Matsubara modes ||. This behavior is confirmed by lattice simulations 
for temperatures that are remarkably close to the critical temperature |J . At lower 
temperatures, the role of higher Matsubara modes become more relevant. One 
purpose of this paper, is to investigate the effects of few Matsubara modes (four) 
in random matrix models. In Section 2, we use the one-flavor version of the NJL 
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Lagrangian in Euclidean space to discuss the effects of mixing between various 
Matsubara modes in random matrix models pj. In Section 3, we specialize to the 
case of Wilson fermions and derive Pastur's equation JllJ for the pertinent random 
matrix model for the two cases where the Matsubara modes are either mixed (model 
I), or decoupled (model II). In Section 4, the Dirac spectral distribution for model 
I is shown to follow from a fifth order algebraic equation (Quinto class), while the 
evolution of the end-points follow from a quartic equation. In Section 5, the Dirac 
spectrum associated to model II is discussed using a superposition of a solution to 
a cubic equation (Cardano class). Some further suggestions and conclusions are 
made in Section 6. In Section 7 we give a short overview of the applicability of 
such models to chiral fermions, and particularly we reproduce the valence quark 
condensate calculated on lattice fL2] ]. Finally, in Section 8 we summarize. 

2. Mapping of the NJL Model into RMM 

The relevance of the Matsubara modes in random matrix models along with their 
mixing is best seen using one flavor NJL Lagrangian in Euclidean space, with anti- 
periodic boundary conditions. Specifically 

U = ^(*7 ' d + imty + g — (V^) 2 + (^ f H5^ , (1) 

where g is a fixed coupling constant. Using the standard bosonization prescription 
P^-2ig 2 4,y L , P^~2ig 2 ^ R , (2) 

where P, P^ stand for independent auxiliary fields, we may rewrite (Q) in the form 

C 4 = +W R {il- d)ip L + ^(i7- d)ip R 

+ ^ R i(P + m)-i/)n + 'ipL i ( pJ ' +m)ip L + —2PPi (3) 

in the chiral basis ip = (i/jr^l)- Note that the Minkowski fields follow from the 
Euclidean fields through (iip^,ip) — > Equation (|J) is defined on the strip 

(3 x Vz in Euclidean space, with P(t + (3, x) = P(t, x), and ip(r + (3, x) = —ip(r, x), 
and a 3-momentum cut-off A. We simplify further the model by putting the left 
and right quark fields on a discrete grid spanned by the three-space points x = 
1,2, ... ,N and introducing the rescaled fields q* = \fVziy^ with 



<Kr,2)= E e-^" T </£, (4) 



n— — oo 



where u> n — (2n+\)irT are the Matsubara frequencies (T = 1/(3). In what follows we 
choose P t— independent. This is equivalent to restricting the four-quark interaction 

/ dT ((^V) 2 + (^75^) 2 ) = 4 E S n+k,m+l ^Rn^Rm^lk^Ll (5) 
JO , , 
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to only the terms with n — m and k = I. Other choices are possible. 

In the case where the auxiliary fields P are r- independent, we can bosonize the 
present model in a different but complementary way ]l0| . After integrating (j^) over 
P and Pt, we bosonize the chirality-flipping pairs using 

which is a doubly banded, complex matrix R with dimensions (N x N) <X> (oo x oo). 
The upper indices refer to three-space x, y and the lower indices to frequency space 
n, m. In terms of (||), the NJL Lagrangian ([!]) becomes 

£o = <? f (n74 + tm)q + ATSTr x>n (RRt) + q* a Bjq L + q\B?qR , (7) 

where the trace in ^ is over x and n. The partition function associated to (Q) is 
simply 

Z[T, fx]=[dR e -^s Tlx ,„(RRt) det2>X)nQ (8) 
with the Dirac operator in a random background, 

where 11 = w n l„ ® l x . The determinant in (g) is over chirality (2), space (x), 
and frequency space (n). This is an example of a chiral random matrix model 
Ul. The structure of (||) generalizes to several flavors as discussed in [||. The 
spectral distribution associated to (^|) exhibits manifest chirality, and is suitable for 
describing Kogut-Susskind fermions on the lattice, as we will discuss below. 



3. Massive Quark Spectra 

The case of Wilson fermions is in a sense different. Indeed, on the lattice, Wilson 
fermions are not manifestly chirally symmetric due to the presence of the r-terms 
full . To compare with lattice Wilson spectra, we consider instead the hermitean 
operator Q = 75 m), where Q is the Dirac Hamiltonian in five-dimensions 
and a single flavor. Here, the Dirac matrix 75 plays the role of the /3-Dirac matrix 
in (3+1) dimensions. A number of simulations, using the analog of Q with Wilson 
fermions on the lattice have been carried out recently by Kalkreutcr [^5| for two- 
color and two-flavor QCD. His results are in qualitative agreement with random 
matrix theory|||4|,[l(|. In 0] the low energy QCD spectra was calculated taking 
into account the lowest Matsubara pair only. The lattice calculation [jl5| however 
shows the appearance of at least one more scale. One of the possible sources of 
another scale is the presence of the next Matsubara pair. 

In terms of chiral random matrix models, the pertinent ensemble for investi- 
gating Wilson fermions for QCD is the Gaussian unitary ensemble for three and 
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more colors and is the Gaussian orthogonal ensemble(l(| for the case of two colors, 
provided that the quarks are in the fundamental representation of color. With this 
in mind, the pertinent random matrix model for spatially constant Wilson modes 
of one-flavor on a cylinder is 

mf d T 



where d T is the r-derivative on a cylinder of length (3 = 1/T, and R a random 
matrix. Equation (|l0| ) is defined with anti-periodic boundary conditions (Q). 

The operator (^0^ , rewritten in the basis (Q) is a the sum of a deterministic and 
a random piece, 

^=(-7n :"WrV (id 



Due to the presence of Wilson r-terms, the Dirac operator is only hermitean, the 
chiral structure is lost. In contrast to Qs in (^), Qw in (0) is not block-off-diagonal. 
The resolvent for the operator (O) is defined as 



= TT^T (Tr-4— ) , ,12) 



2NN* \ z-Q 



w 

where TV* is the number of pairs of Matsubara frequencies retained. The averaging 
depends on the r dependence of the random matrix R, as discussed above in the 
context of the NJL model. Two cases will be considered in this paper : (I) The case 
where R is r-dependent, in which case the various Matsubara modes can mix. (II) 
The case where the matrix R is r-independent, with a block diagonal structure 

R = ®R„ (13) 

in which case, the various Matsubara frequencies are decoupled. 

Throughout, and for simplicity, we choose Gaussian weights V(a) = \a 2 for the 
random distributions 

(...) = ^ J ... cxp[-2iV5] 2 TrV(R)] dR , (14) 

although certainly higher (polynomial) weights are possible jl7|. Note that since 
the dimension of the Gaussian matrix appears explicitly in the random weight, the 
widths of the Gaussians for the mixed and non-mixed case differ and are related by 
the factor ^/N*. In Sections 4 and 5 we use Eq = 1 and T,q = 1 for the Quinto 
and Cardano classes, respectively. In Section 6, in comparing the two classes we set 
2Sq = Y? c = 1. 

The eigenvalue distribution of the Dirac operator ( pi] ) is related to the discon- 
tinuity of G(z) through the real axis 

v{\) = — ImG(z = X + iO). (15) 

7T 
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Using either the law of addition of random matrices [ [L8[ or diagrammatic techniques 
p9| , it follows that the resolvent G(z) is given by the solution of the algebraic 
equation (Pastur equation) 1 1 1 . When the Matsubara frequencies are mixed (model 
I), the corresponding Pastur equation reads 

<*')-wX7=avrZ' (16) 

where are the eigenvalues of the deterministic part of the Dirac operator (^). The 
equation is a polynomial of 2A» + 1 order in G{z). When the Matsubara frequencies 
are decoupled (model II), the corresponding Pastur equation is 

* n—l 

where each of the G n (z) satisfies a cubic (Cardano class) equation 

Gn{z) = 2 {z - G n {z) -M n + z- G n (z) + M„) ' (18) 

We note that M n (T) = \J m 2 + uj 2 , in overall agreement with dimensional reduction 
arguments from high temperature QCD 



4. Coupled Wilson Fermions 



For Wilson fermions, we consider ( |ll| ) in the case of two lowest Matsubara frequen- 
cies. All our considerations, will be made for a single flavor. For the coupled case 
(model I), the Dirac operator reads explicitly 



(19) 



The resolvent of the deterministic part for the two lowest pairs of frequencies (u;_ T 
— ui n -i) reads 





f m 
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with 



Ml = m 2 



((2n + 1)ttT) 5 



(20) 



(21) 



Here n = 0,1 and a (hereafter equal to 1/2) refers to the relative weight of the 
modes. Adding the corresponding Blue's functions of the deterministic and the 
random part one arrives at a fifth order algebraic equation p0| 
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0=0.5 




Fig. 1. Spectral function for Quinto class for the lowest two Matsubara frequencies 
with small mass m = 0.1 in units of Sq = 1 



G 5 + a 4 G 4 + azG 3 + a 2 G 2 + ai G + a = (22) 
with the following coefficients 

0,4 = — 4z 

a 3 = 1 - (M 2 + M 2 ) + 6z 2 

a 2 = -4z 3 -3z + 2z(Af 2 + A/ 1 2 ) (23) 
ai = a(M 2 - Ml) - Af 2 (l - M\) + z 2 [3 - (M 2 + M 2 )] + z 4 
a Q = z [Mq — a(M 2 — M 2 ) — z 2 ] . 

The spectral function generated by the solutions to ( p^) through (^5|) allow for 
a rich phase structure (Quinto class). In general, four possible phases Pi, P 2 , P3, P4 
defined by the number of allowed disconnected arcs or supports of the eigenvalue 
distribution, are possible. In Fig. |l| we show the distribution of eigenvalues as a 
function of temperature and for a single light quark flavor of mass m = 0.1 in units 
of 1 /Yiq . At zero temperature, the spectral function is peaked around zero virtuality, 
with a nonzero condensate. The system is in the Pi phase, and the distribution 
is Wigner's semicircle. As the temperature increases, the condensate decreases, 
followed by the decoupling of the Matsubara modes. The heaviest, the first. The 
structural changes are : P\ — > P 2 Pi or P\ — > P3 — > P4 or Pi — > P4, depending 
on whether the restoration of chiral symmetry precedes, follows or parallels the 
splitting of the frequencies, respectively. Figure || shows the behavior of the spectral 
distribution for a heavier quark flavor. In this case, the quarks are always localized 
with no condensate whatever the temperature. An increase in the temperature 
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Fig. 2. Spectral function for Quinto class for the lowest two Matsubara frequencies 
with large mass to = 2 in units of Eq = 1 

causes only a decoupling of the two Matsubara modes, that is a structural phase 
change from P-z — ► P4. 

Since the general solutions of a fifth order algebraic equation do not exist in 
analytical form, the results shown in Figs. |l]j2|were obtained numerically. Out of the 
five possible solutions, four are ruled out by enforcing normalizability and positivity 
of the spectral distribution. To check that the presented solution is indeed correct 
and unique, we have constructed analytically the equations for the evolution of 
the endpoints. These equations are best amenable in the form of Blue's functions, 
i.e. the functional inverses of the resolvents: G[B(z)] — z. Due to the cuts of the 
Green's functions, the evolution equations stem from the condition B'(z)\ ep = 0, 
where prime denotes differentiation with respect to z and the derivative is taken 
at the endpoint (e.p.). Pastur equations rewritten in terms of the Blue's functions 
are always one degree lower than those for the Green's functions Jig] , therefore 
the evolution equations for the endpoints in our case are solution to a forth order 
(Ferrari) equation. 

In the case considered here, the endpoints ±A; (A\ > A2 > A3 > A4) of the 
four arcs [— A\,— A2], [— A3,— A4], [^4,^3], [^2,Ai] are explicitly given by 




(24) 



F(a, M 2 , Ml) = a(M 2 - M 2 ) - M 2 + M 2 M 2 
where Xi are real positive roots of the Ferrari equation 
X 4 - b 3 X 3 - b 2 X 2 - biX - b = 



(25) 



with coefficients 



1 + 2(Af 2 + Ml) 
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b 2 = 3a(Af 2 - Ml) - (M 2 + Ml) 2 - 2Af 2 (l + M 2 ) 

61 = -a(M ( ^-Af 1 4 ) + 2Af 2 M 1 2 (M 2 + Af 1 2 + l) + M 4 

60 = -aM$Ml{M$-Ml)+M$Ml(l-Ml). (26) 

The critical points are determined from the condition that the real positive points 
At, A3, A2 melt or vanish. In Fig. |^ we show the evolution of the endpoints for 
the numerically generated solutions above (dashed curve) versus the analytically 
generated solutions (solid curve) from (|^). Numerically, the endpoint curve was 
generated using the limit of the norm |ImG(z) — e| — > 0, with e = 0.01 and G(z) 
the normalizable solution to ((2^). The agreement of the two curves, confirms the 
uniqness of the spectral distribution discussed above. 



5. Decoupled Wilson Fermions 

The case of Wilson fermions with two lowest but decoupled Matsubara modes 
(model II), follows from through 

G(z) = ^G Q (z) + ^^-G 1 (z) (27) 
with each G n (z) (n — 0, 1) satisfying the cubic Pastur equation 

Gn{z)= {z-GnW-M*- ^ 

The evolution equations for the endpoints in this case come from the overlap of 
individual supports for each frequency. Since for each frequency the equation for 
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Fig. 4. Spectral function for model II and a light mass m = 0.1 




Fig. 5. Spectral functions for model II and a heavy mass m = 2 



the resolvent is cubic, the corresponding equation for the endpoints is quadratic, 
with the following positive solutions 



1 3 ± y/T+ 8Af2 



V2 1 ± y/l + 8M2 v 

for a — 1/2. Figures |^ and || show the spectral functions versus temperature, for a 
light and heavy quark flavor in units of 1/Sc- 



6. Comparison of the Models 



Despite the overall structural similarities between the spectral functions discussed 
for model I and II, there arc a number of quantitative differences upon closer look. 
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These differences may be helpful for distinguishing these two scenarios in a lattice 
simulation for instance. 

The evolution of the endpoints of the spectra as obtained from model I and II 
are different. Figure ^| shows the behavior of the endpoints versus temperature for 
I (solid line) and II (dashed line) for the light quark flavor case m = 0.1 in units 
where 2Eq = Yjq = 1. The solid curves follow from (^), while the dashed curves 
follow from (pi|). At high temperature (T > 3T C ), the evolution equation becomes 
degenerate. This is expected following the thermal localization. Around the critical 
temperature, however, the evolution of the endpoints are qualitatively different. 

The critical masses and temperatures, are of course different for the two cases 
considered. For model I and at zero virtuality (A = 0) the Quinto equation reduces 
to a second order equation, with the solution 



, E^(M n 2 + M?) - 1 1 i 

G 2 = Q[ ° 2 lJ ±yi + ^ Q (M$-Mty (30) 

for a — 1/2. The support for the spectral density vanishes when G 2 becomes real 
and positive, that is for 

Ml + Ml - 2£ 2 jM 2 M 1 2 = 0. (31) 
This yields the critical temperature 

4^ 2 E 2 T 2 = ^(l-2E 2 m 2 ) + 



mi - 2E 2 m 2 ) 2 - f E 2 m 2 (S 2 m 2 - 1). (32) 



For zero quark mass ( |32| ) gives 7r 2 T c 2 = 5/9 Eg. The critical mass, at which the 
condensate vanishes at zero temperature is m c — 1/Eq. Figure ^ shows the behavior 
of the critical temperature versus m for model I or Quinto class (solid curve). The 
dashed curve is the expected result from model II, with 

7r 2 T 2 = 1 - m 2 . (33) 

For the massless case the critical temperature is in this case ^ , while at zero temper- 
ature the critical mass is 1. The comparison between the two models, was carried 
for 2£q = Eg, = 1. Note that the ratio of the critical parameters depend on the 
mixing scenario 



nT c { 1 model 1 



I model II . ^ 34 ^ 



The spectral distributions around the critical temperature are also qualitatively 
different for the two models considered, as shown in Fig. ||. The spectral function 
in model I following from the Quinto equation (solid curve) is smooth throughout, 
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Fig. 6. Evolution of the endpoints of the spectral distribution for model I (solid 
curve) and model II (dashed curve) with a light quark massm = 0.1 



while the one for model II following from the superposition of two Cardano solutions 
(dashed curve) has sharp edges (discontinuous derivatives). The wiggly behavior 
seen in both spectra is reminiscent of a similar feature seen in Kalkreuter lattice 
spectra with Wilson fermions |{l5[ . This is suggestive of some additional (finite size) 
scale in the latter. This point is worth investigating. 

Finally, in both models the critical exponents are mean- field 1,0,0]. This 
is expected in large N and for Gaussian randomness for all polynomial classes 
(Cardano, Ferrari, Quinto and higher)[| 

7. Staggered Quark Spectra 

For chiral fermions (0) with decoupled Matsubara frequencies, consider the case 
with the dimensional reduction argument Q 

r^II r^n ( ±M n (T) \ ( R \ /orN 

QF = ®±,n Q^n = <8 ± ,» [ ±Mn{T) J + U f J (35) 

with M n (T) = a/ m 2 s + The spectral distribution associated to (|3^) follows from 
the combination of discontinuities of the pertinent solution to a cubic (Cardano) 
equation, Eqs. ( pT|JT^ ). The case considered here at high enough temperature is 
equivalent to the one of Wilson fermions with decoupled Matsubara modes. 

Following Chandrasekharan and Christ ]l^| , let us define the valence quark 
condensate, 

2liir 2mc ^° dA ^' (36) 

1 Quadratic equations do not allow for a phase transition. The system is always in the Pi phase 
(Wigner semicircle). 
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nT=0.9 





Fig. 8. Comparison of the spectral func- 
Fig. 7. Critical temperature as a function tions for model I(solid line), and model II 
of the current mass m for model I (solid (dashed line) with m = 0.1 at nT = 0.9. 
curve) and model II (dashed curve). The The units are fixed by = 2Sq for corn- 
units are fixed by = 2£q = 1 parison 



where ttiq is the valence quark mass, and v{\) the spectral density associated to 
( |35| ) for massless sea quarks, that is m s — 0. Although we have discussed a single 
quark flavor, the analysis carries through unchanged in the massless case. 

In Fig. we show the behavior of the valence quark condensate as obtained by 
Chandrasekharan and Christ jl2], using two-flavor QCD with staggered fermions 
and m s a — 0.01 (or m s — 6 MeV in physical units). In Fig. ^ we display the 
results following from ([56]) with m s = 0, with the two lowest Matsubara pairs (right 
figure) and for comparison the same result but for one Matsubara mode (left figure) 
as discussed in The results are in overall agreement with each other, showing 
mild sensitivity to m s and the addition of further Matsubara modes. The largest 
sensitivity appears around large values of in units of 

Model I for the mixing of the two lowest Matsubara modes is not trivial in the 
case of chiral fermions. Since for the decoupled case model II was equivalent to the 
one of Wilson fermions, we used as model I for chiral fermions the spectral density 
obtained from ( |l9| ) with the proper scaling 2Sq = S^. Comparing the large mass 
part of the spectra we may conclude that the model with the coupling between the 
Matsubara modes agrees better with the lattice calculation (Fig. |9[). 

8. Conclusion 

We have investigated the role of few Matsubara modes on random matrix models 
as inspired from QCD spin and flavor symmetry. After showing how the mixing be- 
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Fig. 9. Semi-quenched condensate versus the valence quark mass for two flavor 
QCD Q with a sea mass m s a = 0.01 
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Fig. 10. Semi-quenched condensate versus the valence quark mass tuq obtained 
from model II (^6|) for the lowest Matsubara modes (left) and the one corresponding 
to the coupling (model I) (right) with zero sea mass 
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tween the various Matsubara modes may affect, the generic structure of the random 
matrix formulation at finite temperature, we have considered two complementary 
models. One with mixing (I) and the other without (II). The spectral function as- 
sociated to the former was shown to follow from a quintic equation (Quinto class), 
while the one for the latter was shown to follow from the superposition of two 
solutions to quartic equations (Cardano class). 

For a light quark mass, the spectral distribution at low temperature has a sup- 
port at zero virtuality. Randomness prevails over temperature. As the temperature 
is increased, the thermal effects cause the spectral distribution to vanish. Ther- 
mal localizations are found to follow the thermal wavelengths. For a heavy quark 
mass, the spectral distribution shows no support at zero virtuality even at zero 
temperature owing to the localization over the Compton wavelength. With increas- 
ing temperature, the thermal localization takes over, with a segregation along the 
thermal frequencies. 

While schematic, the present analysis shows that around the critical tempera- 
ture the effects of the higher Matsubara modes cause quantitative changes not only 
in the critical parameters and their ratios, but also in the endpoints evolution and 
the spectral shapes. These observations are of relevance for lattice simulations of 
Dirac spectra for both Kogut-Susskind and Wilson fermions. 
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